Es werden zwei Experimente durchgeführt:
Nullpunktsbestimmung von Linearpolarisatoren:
Aufbau Nullpunktsbestimmung
Durchführung:
Transmissionsverhalten von Linearpolarisatoren (Variante eine Fiberbench):
Aufbau Transmissionsverhalten von Linearpolarisatoren (A: Egtl. Messung, B: Normierung)
Durchführung:
Transmissionsverhalten von Linearpolarisatoren (Variante zwei Fiberbenches):
Aufbau Transmissionsverhalten von Linearpolarisatoren (A: Egtl. Messung, B: Normierung)
Durchführung:
Nullpunktbstimmung
#
# LINEARPOLARISATOR P1
#
# Extract Data from eLabFTW
P1.maxmin.metadata <- GET.elabftw.byselector(25, node.selector = "#meta-data")[[1]]
P1.maxmin.data <- GET.elabftw.byselector(25, header = T)[[1]]
# Normalise data
P1.maxmin.data$Y1 <- P1.maxmin.data$Y1 / P1.maxmin.metadata[3,2] *100
#
# LINEARPOLARISATOR P2
#
# Extract Data from eLabFTW
P2.maxmin.metadata <- GET.elabftw.byselector(26, node.selector = "#meta-data")[[1]]
P2.maxmin.data <- GET.elabftw.byselector(26, header = T)[[1]]
# Normalise data
P2.maxmin.data$Y1 <- P2.maxmin.data$Y1 / P2.maxmin.metadata[3,2] *100
Winkelabhängige Transmission
#
# LINEARPOLARISATOR P1
#
# Extract data from eLabFTW
P1.transmission.data <- GET.elabftw.byselector(29, header = T)[[1]]
# Normalise data
P1.transmission.data$Y4 <- P1.transmission.data$Y4/P1.transmission.data$Y2 *100
#
# LINEARPOLARISATOR P2
#
# Extract data from eLabFTW
P2.transmission.data <- GET.elabftw.byselector(30, header = T)[[1]]
# Normalise data
P2.transmission.data$Y4 <- P2.transmission.data$Y4/P2.transmission.data$Y2 *100
#
# LINEARPOLARISATOR P3
#
# Extract data from elabFTW
P3.absorbance <- GET.elabftw.bycaption(78, header=T, outputHTTP=T) %>%
parseTable.elabftw(., func=function(x) qmean(x[,4], 0.8, na.rm=T, inf.rm=T),
header=T, skip=14, sep=";") %>%
.[[1]]
colnames(P3.absorbance) <- c("P3", "P4", "background", "measured")
# Normalise data
P3.absorbance$transmittance <- P3.absorbance$measured / P3.absorbance$background
# Plot Zero Point
ggplot(data = P1.maxmin.data, mapping = aes(x = X, y = Y1)) +
geom_point(size = 0.5) +
theme_minimal() +
labs(title = "Die Durchlässigkeit des Linearpolarisators P1 in Abhänigkeit des Rotationswinkels",
x = "Winkel / °",
y = "Anteil der Strahlung, die den Polarisator passiert / %")
# Plot Transmissionsverhalten
ggplot(data = P1.transmission.data, mapping = aes(x = Y3, y = Y4)) +
geom_point() +
theme_minimal() +
labs(title = "Winkelabhängige Transmission von Linearpolarisator P1",
x = "Rotationswinkel / °",
y = "Anteil der Strahlung, die den Polarisator passieren / %")
Fazit:
# Plot Zero Point
ggplot(data = P2.maxmin.data, mapping = aes(x = X, y = Y1)) +
geom_point(size = 0.5) +
theme_minimal() +
labs(title = "Die Durchlässigkeit des Linearpolarisators P2 in Abhänigkeit des Rotationswinkels",
x = "Winkel / °",
y = "Anteil der Strahlung, die den Polarisator passiert / %")
# Plot Transmissionsverhalten
ggplot(data = P2.transmission.data, mapping = aes(x = Y3, y = Y4)) +
geom_point() +
theme_minimal() +
labs(title = "Winkelabhängige Transmission von Linearpolarisator P2",
x = "Rotationswinkel / °",
y = "Anteil der Strahlung, die den Polarisator passieren / %")
Fazit:
# Plot transmission behaviour
# Does the transmittance of the linear polariser change when rotating the polariser and the laser in the same manner?
ggplot(data = P3.absorbance,
mapping = aes(x = P3, y = transmittance*100) ) +
geom_bar(stat="identity") +
theme_classic() +
theme( axis.text = element_text(size=12),
axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1),
panel.grid.major.y = element_line("black", size = 0.1),
panel.grid.minor.y = element_line("grey", size = 0.5) ) +
scale_y_continuous(expand = c(0,0)) +
scale_x_continuous(breaks = P3.absorbance$P3) +
labs(title = expression(bold("The maximal transmittance of the linear polariser P3")),
x = expression(bold("rotation of P3 or the lasers plane of polarisation / °")),
y = expression(bold("transmittance / %"))
)
Fazit:
Es werden zwei Experimente durchgeführt:
Nullpunktbestimmung von Wellenplatten:
Aufbau Nullpunktbestimmung von Wellenplatten (A: Egtl. Messung, B: Normierung)
Durchführung:
Winkelabhängiges Transmissionsverhalten von Wellenplatten:
Aufbau Transmissionsverhalten von Wellenplatten (A: Egtl. Messung, B: Normierung)
Durchführung:
Nullpunkt/Transmissionsverhalten:
#
# Wave Plate W1
#
# Extract data from eLabFTW
W1.data <- GET.elabftw.byselector(27, header = T)[[1]]
# Replace NA values by using the previous value in the vector
fillVector <- function(vector) {
for (index in seq_along(vector)) {
if (is.na(vector[index])) { vector[index] <- vector[index-1] }
}
return(vector)
}
W1.data$Y2 <- fillVector(W1.data$Y2)
W1.data$Y4 <- fillVector(W1.data$Y4)
# Normalise
W1.data$Y1 <- (W1.data$Y1/W1.data$Y2) *100
W1.data$Y3 <- (W1.data$Y3/W1.data$Y4) *100
#
# Wave Plate W2
#
# Extract data from eLabFTW
W2.data <- GET.elabftw.byselector(28, header = T)[[1]]
# Replace NA values by using the previous value in the vector
fillVector <- function(vector) {
for (index in seq_along(vector)) {
if (is.na(vector[index])) { vector[index] <- vector[index-1] }
}
return(vector)
}
W2.data$Y2 <- fillVector(W2.data$Y2)
W2.data$Y4 <- fillVector(W2.data$Y4)
# Normalise
W2.data$Y1 <- (W2.data$Y1/W2.data$Y2) *100
W2.data$Y3 <- (W2.data$Y3/W2.data$Y4) *100
# Plot Transmissionverhalten
ggplot(data = W1.data, mapping = aes(x = X, y=Y3)) +
geom_point() +
theme_minimal() +
labs(title = "Durchlässigkeit der Wellenplatte W1 in Abhängigkeit des Rotationswinkels",
x = "Rotationswinkel / °",
y = "Anteil der Strahlung, die die Wellenplatte passiert / %")
# Plot Nullpunktsbestimmung
ggplot(data = W1.data, mapping = aes(x=X, y=Y1)) +
geom_point(size=0.5) +
theme_minimal() +
labs(title = "Nullpunktsbestimmung der Wellenplatte W1",
x = "Rotationswinkel / °",
y = "Gesteigerte Transmission durch die Wellenplatte / %")
Fazit:
# Plot Transmissionsverhalten
ggplot(data = W2.data, mapping = aes(x = X, y=Y3)) +
geom_point() +
theme_minimal() +
labs(title = "Durchlässigkeit der Wellenplatte W2 in Abhängigkeit des Rotationswinkels",
x = "Rotationswinkel / °",
y = "Anteil der Strahlung, die die Wellenplatte passiert / %")
# Plot Nullpunktbestimmung
ggplot(data = W2.data, mapping = aes(x=X, y=Y1)) +
geom_point(size=0.5) +
theme_minimal() +
labs(title = "Nullpunktsbestimmung der Wellenplatte W2",
x = "Rotationswinkel / °",
y = "Gesteigerte Transmission durch die Wellenplatte / %")
Es werden drei Experimente durchgeführt:
Messung der Stokesvektoren:
Aufbau Messung Stokesvektoren (A: Stokes hinter der Faser, B: Normierung hinter der Faser, C: Stokes vor der Faser, D: Normierung vor der Faser)
Duchführung:
Rotation der Polarisationsebene:
Aufbau Rotationsverhalten optischer Fasern (A: hinter der Faser, B: vor der Faser)
Duchführung:
Definition von Funktionen, die die Messdaten auswerten:
#
# CALCULATE STOKES VECTORS AND THEIR PROPERTIES
#
# Normalise stokes vector and compute polarisation ratio and such shit
# ASSUMPTION: S3 = 0
process.stokesVec <- function(stokes) {
# Compute properties of stokes vectors and normalise -> polarisation ratio, ...
stokes <- lapply(stokes, function(table) {
# Normalise stokes vectors
table[,c("S0", "S1", "S2")] <- table[,c("S0", "S1", "S2")] / table$S0
# Polarisation ratio
table$polarisation <- sqrt(table$S1^2 + table$S2^2) / table$S0
# Polar stokes angle
table$sigma <- better.acos( table$polarisation*table$S0, table$S1, table$S2)
# Polar electrical field coordinate
# !!! Keep in mind: This may be bullshit, if the polarisation ratio is smaller than one!
# I'm to 90% sure that this conversion is also valid for partially polarised light
table$epsilon <- table$sigma / 2
# Return result
return(table)
})
# CALCULATE CHANGE OF THE STOKES VECTORS PROPERTIES
# Change in epsilon, change in polarisation, change in laser intensity
stokes[["change"]] <- data.frame("W" = stokes$PRE$W,
# How much gets the plane of polarisation rotated?
# Calculate the difference as smallest signed distance between two modular numbers
# !!! Keep in mind: This is bullshit, if the polarisation ratio is smaller than one!
"mod.change.in.epsilon" = better.subtraction(stokes$POST$epsilon - stokes$PRE$epsilon, base= pi),
"mod.change.in.sigma" = better.subtraction(stokes$POST$sigma - stokes$PRE$sigma , base=2*pi),
# How much does the polarisation ratio change?
"change.in.polarisation" = stokes$POST$polarisation / stokes$PRE$polarisation - 1,
# How much does the intensity of the light change?
"change.in.intensity" = stokes$POST$I / stokes$PRE$I
)
return(stokes)
}
#
# DO THE STATISTICS
#
# error.stokes contains the data for several identical measurements
# following function will therefore compute statistical properties
# like sd or mean for every column of the tables
# !!! Keep in mind: The statistics for the angles sigma and epsilon are not done with modular functions
# The modular numbers sigma and epsilon will therefore have the wrong mean, variance, sd, ...
# In this case only the modular calculated difference can be trusted.
do.statistics <- function(error.stokes) {
lapply(error.stokes, function(table) {
stats.table <- data.frame( var = sapply(table, var ),
sd = sapply(table, sd ),
mean = sapply(table, mean)
)
return(stats.table)
}) %>% return
# TODO: t-Test, Kruskal-Wallis-Test
}
#
# PLOT THAT SHIT
#
# How does the POLARSATION RATIO change RELATIVE to the initial polarisation ratio?
plot.polarisation.change <- function(data,
title = expression(bold("The Depolarising Behaviour Of <YOUR OPTICAL FIBER>"))
) {
# EXPECTED PARAMETERS:
# data : processedStokesExperiments (output of process.stokesVec)
# title : expression(bold("The Depolarising Behaviour Of <YOUR OPTICAL FIBER>"))
# Create the plot
ggplot(data = data$change,
mapping = aes(x = as.factor(W), y = change.in.polarisation*100) ) +
geom_bar(stat="identity") +
theme_classic() +
theme( axis.text = element_text(size=12),
panel.grid.major.y = element_line("black", size = 0.1),
panel.grid.minor.y = element_line("grey", size = 0.5) ) +
labs(title = title,
x = expression(bold("the have-waveplates angle of rotation "*omega*" / °")),
y = expression(bold("the relative change in the ratio of polarisation "*Delta*Pi*" / %")) )
}
# COMPARING POLARISATION RATIOS before and after interacting with the fiber
plot.polarisation <- function(data,
statistics,
title = expression(bold("The Effect Of An <YOUR OPTICAL FIBER> On The Polarisation Ratio "*Pi))
) {
# PARAMETERS
# data : processedStokesExperiments (output of process.stokesVec)
# statistics : processsedErrorMeasurementExperiment (output of process.stokesVec)
# title : expression(bold("The Effect Of An <YOUR OPTICAL FIBER> On The Polarisation Ratio "*Pi))
ggplot(data = data.frame(W = c(data$POST$W, data$PRE$W),
polarisation = c(data$POST$polarisation, data$PRE$polarisation),
group = c( rep("B_POST", length(data$POST$W)), rep("A_PRE", length(data$PRE$W)) )
),
mapping=aes(x=as.factor(W), y=polarisation*100, fill=group)) +
geom_bar(stat="identity", position = "dodge") +
theme_classic() +
scale_y_continuous(breaks = seq(from=0, to=110, by=10),
expand=c(0,0)) +
theme(axis.text = element_text(size=12),
panel.grid.major.y = element_line("black", size = 0.1),
panel.grid.minor.y = element_line("grey", size = 0.5) ) +
labs(title = title,
x = expression(bold("the have-waveplates angle of rotation "*omega*" / °")),
y = expression(bold("the polarisation ratio "*Pi*" / %")),
fill = "" ) +
scale_fill_discrete( labels=c( expression(bold("before")), expression(bold("after")) ) ) +
geom_errorbar(data=data.frame(W = c(data$PRE$W, data$POST$W) %>% as.factor,
upper = c(data$PRE$polarisation+statistics$POST["polarisation","sd"]*3, data$POST$polarisation+statistics$POST["polarisation","sd"]*3)*100,
lower = c(data$PRE$polarisation-statistics$POST["polarisation","sd"]*3, data$POST$polarisation-statistics$POST["polarisation","sd"]*3)*100,
group = c( rep("A_PRE", length(data$POST$W)), rep("B_POST", length(data$PRE$W)) ),
polarisation = c(data$PRE$polarisation, data$POST$polarisation) ),
mapping = aes(x=W, ymin = lower, ymax=upper, group=group),
position = "dodge"
)
}
# CHANGE in LASER POWER due to optical fiber
plot.intensity.change <- function(data,
title = expression(bold("The Transmittance Of <YOUR OPTICAL FIBER>"))
) {
# EXPECTED PARAMETERS:
# data : processedStokesExperiments (output of process.stokesVec)
# title : expression(bold("The Depolarising Behaviour Of <YOUR OPTICAL FIBER>"))
ggplot( data = data$change,
mapping = aes(x = as.factor(W), y = change.in.intensity*100) ) +
geom_bar(stat="identity") +
theme_classic() +
theme(axis.text = element_text(size=12),
panel.grid.major.y = element_line("black", size = 0.1),
panel.grid.minor.y = element_line("grey", size = 0.5) ) +
labs(title = title,
x = expression(bold("the have-waveplates angle of rotation "*omega*" / °")),
y = expression(bold("the transmitted part of the laser P"[trans]*" / %")) )
}
# COMPARING LASER POWER before and after interacting with the fiber
plot.intensity <- function(data,
title = expression(bold("The Effect Of <YOUR OPTICAL FIBER> On The Lasers Power "*P))
) {
# EXPECTED PARAMETERS:
# data : processedStokesExperiments (output of process.stokesVec)
# title : expression(bold("The Depolarising Behaviour Of <YOUR OPTICAL FIBER>"))
ggplot(data = data.frame(W = c(data$POST$W, data$PRE$W),
intensity = c(data$POST$I, data$PRE$I) / data$PRE$I,
group = c( rep("B_POST", length(data$POST$W)), rep("A_PRE", length(data$PRE$W)) )
),
mapping=aes(x=as.factor(W), y=intensity*100, fill=group)) +
geom_bar(stat="identity", position = "dodge") +
theme_classic() +
theme(axis.text = element_text(size=12),
panel.grid.major.y = element_line("black", size = 0.1),
panel.grid.minor.y = element_line("grey", size = 0.5) ) +
labs(title = title,
x = expression(bold("the have-waveplates angle of rotation "*omega*" / °")),
y = expression(bold("the normalised laser power "*P*" / %")),
fill = "" ) +
scale_fill_discrete( labels=c( expression(bold("before")), expression(bold("after")) ) )
}
# How does the fiber influence the PLANE OF POLARISATIONS ORIENTATION?
plot.plane.rotation <- function(data,
title = expression(bold("The Impact Of The <YOUR FIBER> On The Orientation Of The Plane Of Polarisation"))
) {
# EXPECTED PARAMETERS:
# data : elabFTW table of angle dependent rotation behavior of optical fibers
# title : expression(bold("The Impact Of The <YOUR FIBER> On The Orientation Of The Plane Of Polarisation"))
ggplot( data = data[!is.na(data$X),] ) +
geom_abline( mapping = aes(intercept = 0, slope = 2, color="ideal waveplate") ) +
geom_abline( mapping = aes(intercept = 0, slope = -2, color="ideal waveplate") ) +
geom_point( mapping = aes(x = X, y = Y1-Y1[X==0], color = "after") ) +
geom_point( mapping = aes(x = X, y = Y2-Y2[X==0], color = "before") ) +
theme_classic() +
theme(panel.grid.major = element_line("black", size=0.1),
panel.grid.minor = element_line("grey", size=0.1) ) +
scale_x_continuous(breaks = seq(from=-20, to=300, by=20) ) +
scale_y_continuous(breaks = seq(from=-360, to=360, by=90) ) +
labs(title = title,
x = expression(bold("rotation waveplate / °")),
y = expression(bold("rotation linear polariser / °")),
color = "" )
}
Nutzung der eben definierten Funktionen:
# F1 : POLARISATION-MAINTAING-FIBER
# FETCH data from elabftw
# stokes vector
F1.data.elab <- lapply(c(58, 67, 68), function(experimentID) {
GET.elabftw.bycaption(experimentID, header=T, outputHTTP=T) %>%
parseTable.elabftw(., func=function(x) qmean(x[,4], 0.8, na.rm=T, inf.rm=T),
header=T, skip=14, sep=";")
}) %>% better.rbind(., sort.byrow=1)
# Get the measurements for the error estimation
F1.error.elab <- GET.elabftw.bycaption(66, header=T, outputHTTP=T) %>%
parseTable.elabftw(., func=function(x) qmean(x[,4], 0.8, na.rm=T, inf.rm=T),
header=T, skip=14, sep=";")
# COMPUTE stokes vectors and do statistics on the error estimations
F1.data.stokes <- getStokes.from.expData(F1.data.elab) %>% process.stokesVec
F1.error.stats <- getStokes.from.expData(F1.error.elab) %>% process.stokesVec %>% do.statistics
# F2 : SINGLE-MODE-FIBER
# FETCH data from elabftw
# stokes vectors
F2.data.elab <- lapply(c(69, 70), function(experimentID) {
GET.elabftw.bycaption(experimentID, header=T, outputHTTP=T) %>%
parseTable.elabftw(., func=function(x) qmean(x[,4], 0.8, na.rm=T, inf.rm=T),
header=T, skip=14, sep=";")
}) %>% better.rbind(., sort.byrow=1)
# Rotation of plane of polarisation
F2.rotation.elab <- GET.elabftw.bycaption(72, header=T)[[1]]
# Get the measurements for the error estimation
F2.error.elab <- GET.elabftw.bycaption(71, header=T, outputHTTP=T) %>%
parseTable.elabftw(., func=function(x) qmean(x[,4], 0.8, na.rm=T, inf.rm=T),
header=T, skip=14, sep=";")
# COMPUTE stokes vectors and do statistics on the error estimations
F2.data.stokes <- getStokes.from.expData(F2.data.elab) %>% process.stokesVec
F2.error.stats <- getStokes.from.expData(F2.error.elab) %>% process.stokesVec %>% do.statistics
# FIBER F3 : MULTI-MODE-FIBER
# FETCH data from elabftw
# stokes vectors
F3.data.elab <- GET.elabftw.bycaption(74, header=T, outputHTTP=T) %>%
parseTable.elabftw(., func=function(x) qmean(x[,4], 0.8, na.rm=T, inf.rm=T),
header=T, skip=14, sep=";")
# Rotation of plane of polarisation
F3.rotation.elab <- GET.elabftw.bycaption(73, header=T)[[1]]
# Get the measurements for the error estimation
F3.error.elab <- GET.elabftw.bycaption(75, header=T, outputHTTP=T) %>%
parseTable.elabftw(., func=function(x) qmean(x[,4], 0.8, na.rm=T, inf.rm=T),
header=T, skip=14, sep=";")
# COMPUTE stokes vectors and do statistics on the error estimations
F3.data.stokes <- getStokes.from.expData(F3.data.elab) %>% process.stokesVec
F3.error.stats <- getStokes.from.expData(F3.error.elab) %>% process.stokesVec %>% do.statistics
# How does the polarisation ratio change relative to the initial polarisation ratio?
plot.polarisation.change(data = F1.data.stokes,
title = expression(bold("The Depolarising Behaviour Of An Optical PM-Fiber (F1)") )
)
# COMPARING POLARISATION RATIOS before and after interacting with the fiber
plot.polarisation(data = F1.data.stokes,
statistics = F1.error.stats,
title = expression(bold("The Effect Of An Optical PM-Fiber (F1) On The Polarisation Ratio "*Pi))
)
# CHANGE in LASER POWER due to optical fiber
plot.intensity.change(data = F1.data.stokes,
title = expression(bold("The Transmittance Of An Optical PM-Fiber (F1)"))
)
# COMPARING LASER POWER before and after interacting with the fiber
plot.intensity(data = F1.data.stokes,
title = expression(bold("The Effect Of An Optical PM-Fiber (F1) On The Lasers Power "*P))
)
Fazit:
# How does the polarisation ratio change relative to the initial polarisation ratio?
plot.polarisation.change(data = F2.data.stokes,
title = expression(bold("The Depolarising Behaviour Of An Optical SM-Fiber (F2)") )
)
# COMPARING POLARISATION RATIOS before and after interacting with the fiber
plot.polarisation(data = F2.data.stokes,
statistics = F2.error.stats,
title = expression(bold("The Effect Of An Optical SM-Fiber (F2) On The Polarisation Ratio "*Pi))
)
# CHANGE in LASER POWER due to optical fiber
plot.intensity.change(data = F2.data.stokes,
title = expression(bold("The Transmittance Of An Optical SM-Fiber (F2)"))
)
# COMPARING LASER POWER before and after interacting with the fiber
plot.intensity(data = F2.data.stokes,
title = expression(bold("The Effect Of An Optical SM-Fiber (F2) On The Lasers Power "*P))
)
# How does the fiber influence the PLANE OF POLARISATIONS ORIENTATION
plot.plane.rotation(F2.rotation.elab,
title = expression(bold("The Impact Of The Single-Mode Fiber (F2) On The Orientation Of The Plane Of Polarisation"))
)
## Warning: Removed 4 rows containing missing values (geom_point).
Fazit:
# How does the polarisation ratio change relative to the initial polarisation ratio?
plot.polarisation.change(data = F3.data.stokes,
title = expression(bold("The Depolarising Behaviour Of An Optical MM-Fiber (F3)") )
)
# COMPARING POLARISATION RATIOS before and after interacting with the fiber
plot.polarisation(data = F3.data.stokes,
statistics = F3.error.stats,
title = expression(bold("The Effect Of An Optical MM-Fiber (F3) On The Polarisation Ratio "*Pi))
)
# CHANGE in LASER POWER due to optical fiber
plot.intensity.change(data = F3.data.stokes,
title = expression(bold("The Transmittance Of An Optical MM-Fiber (F3)"))
)
# COMPARING LASER POWER before and after interacting with the fiber
plot.intensity(data = F3.data.stokes,
title = expression(bold("The Effect Of An Optical MM-Fiber (F3) On The Lasers Power "*P))
)
# How does the fiber influence the PLANE OF POLARISATIONS ORIENTATION
plot.plane.rotation(F3.rotation.elab,
title = expression(bold("The Impact Of The Multi-Mode Fiber (F3) On The Orientation Of The Plane Of Polarisation"))
)
Fazit:
Es wird ein Experiment mit zwei verschiedenen Aufbauten durchgeführt:
Polarisationsempfindlichkeit:
Aufbau Polarisationsempfindlichkeit des WiTec-Detektors (A: Mit Mikroskop, B: Ohne Mikroskop)
Durchführung:
# CONVERT TIME SERIES OF SPECTRA INTO EASY PLOTABLE DATA.FRAME
# This function turns a data.frame with multiple spectra organised in multiple columns into a data.frames with all spectra staked in the same columns
# Makes plotting with ggplot easier
makeSpectraPlotable <- function(spectra, colorFunc=function(x) return(x)) {
# PARAMETERS
# spectra : the return value of parseTimeSeries.elab, contains a time series with multiple spectra
# colorFunc : a custom function applied to the time variable (in this case the position of the linear polariser) to tweak the color scale
# RETURN VALUE
# A data.frame with all spectra stacked on top of each other. data.frame has the columns: wavenumber, signal, P (linear polarisers rotation), color
# Extract the position of the linear polariser from the column names and repeat the value to match the length of the corresponding spectrum
# Used for colouring and grouping the data correctly when plotting with ggplot
polariser <- lapply(seq_along(spectra[,-1])+1, function(index) {
rep( colnames(spectra)[index], length.out=length(spectra[,index]) )
} ) %>% unlist %>% as.numeric
# Create a dataframe with all spectra stacked on top of each other, instead of every spectrum in an own column
data.frame(wavenumber = spectra$wavenumber,
signal = unlist(spectra[,-1]) %>% unname,
P = polariser,
# Hand the polariser position to a custom function to tweak the color scale
color = colorFunc(polariser)
)
}
#
# PLOTING FUNCTIONS
#
# Plot all measured WHITE LIGHT SPECTRA in one plot and color code them according to the rotation of the linear polariser
# The color should show the absolute DEVIATION of the lasers plane of polarisation FROM THE DETECTORS MOST SENSITIVE AXIS
plot.detector.whitelamp <- function(data,
title = "The Changing Detector Response For Different Linear Polarised White Light <Of Your Equipment>"
) {
# PARAMETERS
# data : plotable time series of spectra. Use the return value of RHotStuff::parseTimeSeries.elab() %>% makeSpectraPlotable()
# title : Some descriptive title
ggplot( data = data,
mapping = aes(x = wavenumber, y = signal, group = P, color = color)
) +
scale_color_gradient(low = "blue",
high = "red",
breaks = seq(from=0, to=90, by=22.5)
) +
theme_hot() +
labs(title = title,
y = "counts",
x = expression(bold("wavenumber / cm"^"-1")),
subtitle = "the color gradient encodes the absolute deviation D of the linear polarisers position \nfrom the detectors most sensitive axis",
color = "D / °") +
geom_line()
}
# Plot WHITE LAMP SPECTRA in one 3d plot as 3D SURFACE (single picture)
plot.detector.allSpectra <- function(data,
title = expression(bold("The White Lamp Raman Spectra For Different Polarised Light")),
color.resolution = 100,
color.ramp = c("blue", "red"),
theta = 270,
phi = 20,
grid.resolution.X = 20,
grid.resolution.Y = 2
) {
# Seperate wavenumber axis, polariser position and spectra
PlotMat <- as.matrix(data[, -1])
wavenumber <- data$wavenumber
polariser <- as.numeric( colnames(PlotMat) )
# Create a grid for plotting
grid <- list(ordinate = wavenumber, abcissa = polariser)
grid.surface <- make.surface.grid(grid)
# Create a 3d plottable surface
surface <- as.surface(grid.surface, PlotMat)
# Create color palette
col.Palette <- colorRampPalette(color.ramp)(color.resolution)
# Calculate Color of the surface according to the z-value of the corresponding point
zfacet <- PlotMat[-1, -1] + PlotMat[-1, -ncol(PlotMat)] + PlotMat[-nrow(PlotMat), -1] + PlotMat[-nrow(PlotMat), -ncol(PlotMat)]
facetcol <- cut(zfacet, color.resolution)
plotCol <- persp(surface, theta=theta, phi=phi)
# Create the plor
plot.surface(surface, type="p", theta=theta, border=NA, phi=phi,
xlab = "wavenumber / cm^-1",
ylab = "wave plate position / °",
zlab = "detector signal / counts",
main = title)
# Add grid lines
# Get the position of the gridlines
select.X <- seq(1,length(grid[[1]]), by=grid.resolution.X)
select.Y <- seq(1,length(grid[[2]]), by=grid.resolution.Y)
xGrid <- grid[[1]][select.X]
yGrid <- grid[[2]][select.Y]
# Draw the gridlines
for(i in select.X) lines(trans3d(x=rep(grid[[1]][i],ncol(PlotMat)),
y=grid[[2]],
z=PlotMat[i,], pmat=plotCol))
for(i in select.Y) lines(trans3d(x=grid[[1]],
y=rep(grid[[2]][i],nrow(PlotMat)),
z=PlotMat[,i], pmat=plotCol))
}
# Fetch experimental data from elabFTW
detector.spectra <- GET.elabftw.bycaption(76, header=T, outputHTTP=T) %>% parseTimeSeries.elab(., header=F, sep="")
#
# PLOT THAT SHIT
#
# Plot the white lamp spectra for the detector without the microscope
plot.detector.whitelamp(data=makeSpectraPlotable(detector.spectra[[2]],
colorFunc=function(polariserRotation) {mod(polariserRotation, 180) %>% `-`(.,90) %>% abs(.)} ),
title="The Changing Detector Response For Different Linear Polarised White Light Of The WiTecs Detector")
# Plot the white lamp spectra for the detector with the microscope
plot.detector.whitelamp(data=makeSpectraPlotable(detector.spectra[[1]],
colorFunc=function(polariserRotation) {mod(polariserRotation, 180) %>% `-`(.,90) %>% abs(.)} ),
title="The Changing Detector Response For Different Linear Polarised White Light Of The WiTecs Detector And Microscope")
# How does the sensitivity of the detector change with the wavenumber (quotient method)
ggplot( data = data.frame(wavenumber = detector.spectra[[2]]$wavenumber,
quotient = apply(detector.spectra[[2]][,-1], 1, function(slice) { max(slice)/min(slice) })
),
mapping = aes(x=wavenumber, y=quotient) ) +
theme_hot() +
labs(title = "Quotient of maximal and minimal detector response",
x = expression(bold("wavenumber / cm"^"-1")),
y = "max/min") +
geom_line()
# Plot the WHITE LAMP SPECTRA in one 3d plot as 3D SURFACE
plot.detector.allSpectra(detector.spectra[[2]][,-c(21:24)], theta=240)
Fazit:
Aufbau Messung von Ramanspektren mit spezifischer Laserpolarisation
Durchführung:
# FETCH DATA
tetra.spectra <- GET.elabftw.bycaption(79, header=T, outputHTTP=T) %>%
parseTimeSeries.elab(., col.spectra=3, sep="") %>% .[[1]]
# PREPROCESS
# Statistical Background Correction
# Peaks no longer available
# tetra.spectra[,-1] <- sapply(tetra.spectra[,-1], function(spec) spec-Peaks::SpectrumBackground(spec))
tetra.spectra[,-1] <- t( as.matrix(tetra.spectra[,-1]) ) %>%
baseline::baseline(., method="fillPeaks", lambda=1, it=10, hwi=50, int=2000) %>%
baseline::getCorrected(.) %>% t(.)
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
# Normalise each spectra by its highest peak
tetra.spectra[,-1] <- sapply(tetra.spectra[,-1], function(spec) spec/max(spec))
# FIND THE LOCATION OF THE PEAKS
# Guess the general area (wavenumbers) where the peak lays in
tetra.peakMargins <- list( c(100, 250),
c(250, 340),
c(340, 580),
c(580, 770),
c(770, 890) )
# Find the maximum in each guessed area
tetra.peakLocations <- sapply(tetra.peakMargins, function(margins) {
# Guess the broad location of the peak
selectGeneralPeakLocation <- which(tetra.spectra$wavenumber > margins[1] & tetra.spectra$wavenumber < margins[2])
# Find the maximum of the guessed area and returns its location
peak.value <- max( tetra.spectra$`0`[selectGeneralPeakLocation] )
peak.location <- tetra.spectra$wavenumber[tetra.spectra$`0` == peak.value]
return(peak.location)
})
# EXTRACT THE PEAK HEIGHTS FOR DIFFERENT LASER POLARISATIONS
# Create empty matrix
tetra.peakChange <- matrix( rep(NA, ncol(tetra.spectra[,-1])*(length(tetra.peakLocations)+1) ), nrow = ncol(tetra.spectra[,-1]) )
# Add column with the wave plate position
tetra.peakChange[,1] <- colnames(tetra.spectra[,-1]) %>% as.numeric
# Loop over all peaks and and add their heights to the matrix
for (index in seq_along(tetra.peakLocations)) {
# Get the wavenumber of the peak
peak <- tetra.peakLocations[index]
# Add the corresponding row of the spectra table to the matrix
tetra.peakChange[,index+1] <- tetra.spectra[ which(tetra.spectra$wavenumber==peak),-1 ] %>% unlist
}
# Add describtive column names
colnames(tetra.peakChange) <- c("waveplate", tetra.peakLocations)
# Translate matrix to data.frame
tetra.peakChange <- as.data.frame(tetra.peakChange)
# COMPUTE DETECTORS SENSITIBLITY FOR LIGHT POLARISED ALONG DIFFERENT OPTICAL AXIS
tetra.sensibility <- sapply(tetra.peakChange[,-1], function(peakheight) max(peakheight)/min(peakheight))
tetra.sensibility <- data.frame(wavenumber = names(tetra.sensibility) %>% as.numeric,
quotient = tetra.sensibility %>% unname)
# HWO ARE THE OPTICAL AXIS OF THE DETECTOR ALIGNED?
# The wave plates position of the maximal detector response
# Loop over all peaks and compare the height between the different spectra
sapply(seq_along(tetra.peakChange[,-1]), function(index) {
# The maximal height of the current peak
peakHeight <- max(tetra.peakChange[,index+1])
# Get the index of the maximum
select <- which(peakHeight==tetra.peakChange[,index+1])
# If there are multiple maxima its propably the peak I normalised with
# Just return NA in this case
if (length(select)>1) return(NA)
# Return the wave plates position
return( tetra.peakChange[ select,1 ] )
})
## [1] 158 158 NA 158 45
# Plot ALL SPECTRA OVERLAYERD as 2d plot
tetra.plot.allSpetra <- ggplot( data = makeSpectraPlotable(tetra.spectra[tetra.spectra$wavenumber>100,-c(21:24)],
colorFunc=function(waveplateRotation)
{ mod(waveplateRotation, 90) %>% `-`(., 45) %>% abs } ),
mapping = aes(x = wavenumber, ymax = signal, ymin=0, group = P, fill = color) ) +
scale_fill_gradientn(colors = c("red", "orange", "green"),
breaks = seq(from=0, to=45, length.out=4) ) +
theme_hot() +
labs(title = "Influence Of Light Polarisation On Raman Spectrum Of Tetrachloromethane",
y = "normalised intensity",
x = expression(bold("wavenumber / cm"^"-1")),
subtitle = "the color gradient encodes the absolute deviation D of the wave plates position \nfrom the detectors least sensitive axis",
fill = "D / °") +
geom_ribbon()
# Plot all wavenumbers
tetra.plot.allSpetra
# Show just the interesting part
tetra.plot.allSpetra + coord_cartesian(xlim = c(150, 850), ylim = c(0,0.55))
# Plot the HEIGHT OF PEAKS against the wave plates position
plot(x=tetra.peakChange$waveplate, y=tetra.peakChange[,2], type="n", ylim=c(0,1),
main = "Peakheight Change Due To Polarisation Change",
xlab = "waveplate rotation / °",
ylab = "normalised peak height")
for (index in seq_along(tetra.peakChange[,-1])) lines(x=tetra.peakChange$waveplate, y=tetra.peakChange[,index+1], col=index)
# Plot the quotient of the DETECTOR RESPONSE along its optical axis
plot(tetra.sensibility, type="h", lwd=2,
main = "Quotient of maximal and minimal detector response",
xlab = "max/min",
ylab = expression("wavenumber / cm"^"-1") )
Plot/Diskussion